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Abstract 

The single-step one-shot method has proven to be very efficient for PDE- 
constrained optimization where the partial differential equation (PDE) is 
solved by an iterative fixed point solver. In this approach, the simulation 
and optimization tasks are performed simultaneously in a single iteration. 
If the PDE is unsteady, finding an appropriate fixed point iteration is non¬ 
trivial. In this paper, we provide a framework that makes the single-step 
one-shot method applicable for unsteady PDEs that are solved by classi¬ 
cal time-marching schemes. The one-shot method is applied to an optimal 
control problem with unsteady incompressible Navier-Stokes equations that 
are solved by an industry standard simulation code. With the Van-der-Pol 
oscillator as a generic model problem, the modihed simulation scheme is fur¬ 
ther improved using adaptive time scales. Finally, numerical results for the 
advection-diffusion equation are presented. 

Keywords: simultaneous optimization, one-shot method, PDE-constrained 
optimization, unsteady PDE, adaptive time scale 


1. Simultaneous PDE-constrained optimization 

Many engineering and science problems can be described by partial dif¬ 
ferential equations (PDEs). Instead of solving them analytically, the rapid 
progress in computer technologies made it possible to compute approximate 
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solutions numerically. In Computational Fluid Dynamics (CFD), advanced 
and mature simulation tools have been developed, refined and validated over 
years. Nowadays, high fidelity simulation tools are available that compute 
accurate approximations for a variety of complex fluid flow configurations 

BIZI- 

The simulation task refers to the process of solving the PDF numerically 
for some state variables (as for example the velocity, density or tempera¬ 
ture of a system) while some appropriate data is given (such as geometry, 
material coefficients, boundary or initial conditions). In contrast, the task 
of optimization is to adjust some of this data in such a way, that the state 
variables exhibit a desired behavior determined by an objective function J. 
Let y be the vector of state variables of the system and let u describe the 
data that can be adjusted, the so-called design variables. The optimization 
problem under consideration then reads 

min J{y,u) s.t. c{y,u) = 0 (1) 

where c represents a system of PDFs including boundary and/or initial con¬ 
ditions. Depending on the choice of the state and design variables, the opti¬ 
mization problem can represent an optimal shape design problem, an inverse 
design problem for parameter estimation or an optimal control problem. It 
has a wide range of application scenarios as for example in aerodynamic 
shape design, where one aims to find an airfoil shape that minimizes its drag 
coefficient while the fluid flow satisfies the Navier-Stokes equations, or appli¬ 
cations in other disciplines like geophysics, medical imaging and atmospheric 
science. 

If evaluating the objective function is rather time consuming, the op¬ 
timization problem is typically solved using gradient-based methods Pi]- 
These methods iteratively perform design updates utilizing the sensitivity 
of the objective function to the design variables. If the dimension of the 
design space is rather large, the adjoint method is preferred since the cost 
for computing the sensitivities is then independent of the dimension of the 
design space PEj. In this approach, an adjoint PDF system is derived from 
the Lagrangian function associated with ([^. The solution of the adjoint sys¬ 
tem can be used to determine the desired sensitivities efficiently. Two main 
approaches are distinguished namely the discrete adjoint approach, where a 
discrete adjoint system is derived from the discretized PDFs, and the contin¬ 
uous adjoint approach that derives the adjoint system from the continuous 
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PDEs and discretizes afterwards. The discrete adjoint solution can be gen¬ 
erated automatically using Automatic Differentiation (AD) [7j. 

In most application scenarios it is reasonable to assume, that every design 
variable u uniquely determines a state y{u) that satishes the state equation 
c{y{u),u) = 0. Under this assumption it is a common approach to eliminate 
the state variable from the optimization and focus on the unconstrained 
minimization problem 


min J{y{u),u). (2) 

U 

In this so called reduced space approach a gradient-based optimization strat¬ 
egy can be applied to the reduced objective function ([^ that depends on 
the design u solely, while the PDE-constraint is treated implicitly by recov¬ 
ering c{y{u),u) = 0 after each design change. Methods of this type are also 
referred to as black-box approach or Nested Analysis and Design (NAND) 
[8]. The main drawback of reduced space methods is a direct consequence 
of the implicit treatment of the constraint: After each design change, the 
state variable has to be recomputed such that it satisfies the PDE. This 
means, that a full numerical simulation has to be performed after each de¬ 
sign change. Despite of the rapid growth in computer capacities, simulating 
nonlinear PDE-systems still can take hours or even weeks on state-of-the- 
art supercomputers. This makes the reduced space approach unaffordable in 
many sophisticated application scenarios. 

On the other extreme, so called full space methods are a popular alter¬ 
native |9]. Instead of reducing the optimization space by recovering the 
PDE-solution after each design change, the optimization problem is solved 
in the full space for the (discretized) state and design variables. In this 
approach, the optimality conditions for the constrained minimization prob¬ 
lem ([^ are solved simultaneously for the state, the adjoint and the design 
variable in a SQP-like fashion. Because the simulation is directly integrated 
in the optimization process, these methods are often called Simultaneous 
Analysis and Design (SAND), all-at-once approach, or one-shot approach 
[TUI [TT| [TU| [TUI dH [TUI CS] • It has been observed numerically - at least for 
steady state PDEs - that the full space methods can outperform the re¬ 
duced space methods by about one order of magnitude measured in iteration 
counts and runtime mm- However, the major drawback of many full space 
optimization methods is, that they require the computation of additional Ja- 
cobians and Hessians which are not necessarily part of the PDE simulation 
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tool. In fact, many PDE solvers only approximate Jacobians of the dis¬ 
cretized PDE residuals due to implementation or computational effort which 
makes it necessary to rewrite and enhance the PDE solver for optimization. 

In this paper, we consider a full space optimization method that over¬ 
comes this difficulty. We follow the single-step one-shot approach as proposed 
in Ha [15] that is specially tailored for PDE-constrained optimization where 
it is impossible to form or factor the Jacobian of the constraint. Instead, it 
is assumed that the user is provided with an iterative hxed point algorithm 
that computes a discrete numerical approximation to the PDE solution in 
a black box fashion. In the considered one-shot approach, these iterations 
are enriched by an iteration for the adjoint as well as the design variable, 
so that in each optimization iteration, only one step of the PDE solver and 
the adjoint solver is executed. The iterative adjoint solver can be computed 
efficiently with the reverse mode of AD applied to the PDE solver and eval¬ 
uating the objective function. The design step is based on the gradient of 
the reduced space objective function. Provided that a certain preconditioner 
for the design update is used, the single-step one-shot method is proven to 
converge to an optimal point of the minimization problem [18]. Since the 
iterations of the PDE hxed point solver are used in a black-box manner, the 
optimization method leverages and retains the software investment that has 
been made in developing the PDE solver. Section shortly recalls the main 
aspects of the considered one-shot approach. 

Application of the single-step one-shot method to optimization with steady 
state PDEs is straightforward in terms of the hxed point solver: It is a com¬ 
mon and well established approach for solving steady state PDEs to apply 
the so-called pseudo-time-stepping method. In this approach, the PDE is 
interpreted as a steady state of a dynamical system and solved numerically 
by an explicit (pseudo-)time-stepping method [191 [IB]- Its strong relation to 
general iterative methods made it possible to apply the proposed one-shot 
method to various optimization tasks with steady state PDEs especially in 
computational huid dynamics [m 1201 [I5112I]. 

However, if the PDE is fully unsteady, hnding an appropriate hxed point 
iteration that hts in the single-step one-shot framework raises complexity. 
Existing simulation tools for unsteady PDEs typically apply an implicit time 
marching scheme. The resulting implicit equations are solved one after an¬ 
other forward in time utilizing an iterative solver at each time step as pro¬ 
posed in the famous dual time-stepping approach by Jameson [22] • In order 
to prepare for single-step one-shot optimization, a hxed point iterator for the 
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unsteady PDE is derived from such a method in Section by reducing the 
number of iterations at each time step. Since the iteration steps themselves 
are not changed, the effort that has been spent for developing the PDE-solver 
is preserved within the new scheme. 

In Section the proposed one-shot method is applied to an optimal 
control problem with unsteady incompressible Navier-Stokes equations that 
are solved by an industry standard simulation code. The test case under 
consideration is an active flow control problem of a cylinder in unsteady 
flow. Actuation slits are installed on the cylinder surface in order to reduce 
vorticity in the wake. 

The modihed time marching scheme is further improved using adaptive 
time scales in Section Numerical tests on model problems have shown, 
that the number of iterations needed to solve the unsteady PDE can thereby 
be reduced drastically |23]. The adaptive time scale approach is further 
investigated in Section for the advection-diffusion equation with periodic 
boundary condition. 

2. The single-step one-shot method 

Let c(|/, u) = 0 with c: Y xU ^ Y represent a system of PDEs with state 
vector y and a set of design variables u G U. In a discretized PDE setting, it 
is reasonable to assume that Y and U are hnite dimensional Hilbert spaces 
with dim Y = m and dim U = n which allows us to associate their elements 
with the corresponding coordinate vectors in M™' and respectively, and 
write duals as transposed vectors. For an objective function J\ Y x f/ —)■ M 
we want to solve the PDE-constrained optimization problem 

min J{y,u) s. t. c{y,u)=0. (3) 

y,u 

The proposed one-shot method is tailored for optimization problems where 
the PDE-constraint is solved by an iterative hxed point solver that computes 
a feasible state variable given appropriate design data. We therefore assume, 
that an iteration function H: Y x f/ —)■ T is available which - for any given 
design u* & U - iteratively updates the state variable with 


yk+i = H{yk,u*) 


(4) 
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such that the limit y* = limk^aoVk ^ h" exists and satishes c{y*,u*) = 0. 
Convergence of the above iteration is assured if 


dH 

dy 


< p<l 


(5) 


holds for all points of interest, i.e. the iteration function H is contractive 
with respect to the state variable according to Banach’s fixed point theorem 
j24j . If p is close to 1, the simulation code is a rather slowly converging fixed 
point iteration as it is typically the case for applications in CFD. 

Assuming that the limit y* is a hxed point of the iteration function H if 
and only if c{y*,u*) = 0, we can reformulate the constraint of the optimiza¬ 
tion problem in terms of the hxed point equation y* = H{y*,u*) and focus 
on the following minimization problem 


min J{y,u) s. t. y = H{y,u). (6) 

y,u 

We dehne the associated Lagrangian function 

L{y,y,u) := J{y,u) + iH{y,u) -yfy (7) 


with Lagrange multiplier y G Y* which corresponds to the so-called adjoint 
variable. Any local optimizer of (|^ is a saddle-point of L leading to the 
following necessary optimality conditions, the so-called Karush-Kuhn-Tucker 
(KKT) conditions [3]: 


y = H{y,u) 

state equation 

(8) 

y = Jyiy, + Hy{y, u)'^y 

adjoint equation 

(9) 

0 = Ju{y, -h Huiy, uYy 

design equation 

(10) 


where subscripts denote partial derivatives. The state equation, also referred 
to as the primal equation, can be solved with the hxed point iteration Q 
by construction. The adjoint equation is linear in the adjoint variable y 
and involves the transpose of the Jacobian of the primal hxed point iterator 
with respect to the state. An iteration for solving the adjoint equation can 
be obtained by applying AD to H and evaluating J [7]. This approach 
automatically generates an iteration of the following type 

Tjk+i = Jy{y, uV + Hyiy, ufyk- (11) 
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Since H is contractive, this iteration converges to a solution of (|^ for any 
given design and corresponding state that satisfy the primal equation. The 
design equation refers to stationarity of the objective function with respect 
to design changes. For feasible state and adjoint variables, it corresponds to 
the total derivative of reduced space objective function in (|^ and is often 
called the reduced gradient. The reduced gradient is therefore utilized for 
updating the design variable during the optimization procedure. Again, all 
partial derivatives can be evaluated efficiently with the use of AD. 

Instead of solving the primal equation with Q hrst, then iterating for 
an adjoint solution with ( [IT| ), followed by an update of the design in the 
direction of the reduced gradient, the main idea of the single-step one-shot 
method is to solve the set of KKT-equations simultaneously in one single 
iteration [ElEI]: 


yk+1 


H{yk,Uk) 

Vk+l 

= 

Hyi^Z/k-jUk) Vk 

^fc+i 


'^k (^Jui^yk') ^/c) “1“ ^ 


In this approach, the adjoint iteration and a preconditioned reduced gradient 
step for the design variable are integrated into the primal iteration. In order 
to ensure convergence of the single-step one-shot iteration, the preconditioner 
Bk has to be chosen such that the coupled iteration is contractive. It is 
suggested in [TSl El] to look for descent on the doubly augmented Lagrangian 
function 

L'"{y,y,u) ■= I \\H{y,u) -y\f + ^ \\jy{y,u) + {Hy{y,u) - I^yf 

+ L{y,y,u) (13) 

where weighted residuals of the state and the adjoint equations are added 
to the Lagrangian function with weights a,/9 > 0. It is proven ibidem, that 
a suitable preconditioner approximates the Hessian of the augmented 
Lagrangian with respect to the design u. Numerically, Bk can therefore 
be approximated using secant updates on the gradient of the augmented 
Lagrangian VuL°', as for example applying BFGS-updates in each iteration 
[3]. Computation of the gradient can be automated with the use of AD, 
where forward over reverse differentiation provides routines for computing 
second derivatives Eg. 
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The above single-step one-shot method npdates the primal, the adjoint 
and the design variable simnltaneously in a Jacobi-like fashion. In [26] vari¬ 
ants of the method are snrveyed which update the variables in a Seidel-type 
iteration where the new variables are used immediately. The design space 
preconditioner should then approximate the reduced space Hessian of the 
objective function. Furthermore, multi-step one-shot methods are analysed 
where not only one but several steps of the primal and the adjoint iteration 
are performed before the design is updated. 

The single-step one-shot method has been applied to various optimization 
problems where the underlying PDF is steady [I7ll2Ill2n]. In that case, the 
hxed point iterator H arises naturally in common simulation tools. It can 
be for example one step of an explicit pseudo-time-stepping scheme which 
is typically state of the art for solving steady state PDFs in industrial CFD 
applications DEI Numerical applications have proven, that the cost for an 
one-shot optimization is only a small multiple of the cost of a single simulation 
of the underlying PDF - a property which is called hounded retardation. The 
factor typically varies between 2 and 8. Direct comparison with a classical 
reduced space BFGS optimization showed, that the overall runtime for the 
one-shot optimization is about one order of magnitude lower than for the 
reduced space approach [23 120] . 

3. Single-step one-shot optimization with unsteady PDFs 

With the rapid increase in computational capacities, numerical simula¬ 
tion codes are no longer restricted to steady state solutions but perform 
accurate high hdelity simulations of unsteady turbulent flow. Common sim¬ 
ulation methods discretize the transient term of the PDF in time applying 
a time-stepping method (e.g. Runge-Kutta) while implicit methods are of¬ 
ten preferred due to good stability properties [28]. The resulting implicit 
equations are then solved iteratively as for example in dual time-stepping 
methods with Picard- or Newton-iterations in each time step, often in com¬ 
bination with multigrid methods and Krylov subspace techniques [22l |2^ . 
In this section we present a framework for single-step one-shot optimization 
with unsteady PDFs utilizing these time-marching schemes. 

3.1. Problem statement 

For time-dependent PDFs the state variable varies with time and thus is 
a function y: M. ^ Y. The objective function to be minimized is typically 


given by a time averaged quantity 


J{y,u) 



J{y{t),u) dt 


(14) 


where J: F x —)■ M represents some time-dependent quantity as for exam¬ 
ple drag or lift coefficient of an airfoil in unsteady flow. The optimization 
problem with unsteady PDE-constraints then reads 


min 

y,u 


f j i(|/(t),M)dt 

subject to 

(15) 

dy{t) f, f.-s X 

Vt e [0,T] 

(16) 

1/(0) = y° 


(17) 


where the right hand side f:Y x U Y corresponds to spatial derivative 
operators and boundary terms of the unsteady PDE and y*’ G E is some 
appropriate initial data. 

3.2. Fixed point iteration for unsteady PDE-constraints 

Numerical methods for solving unsteady PDEs discretize the time do¬ 
main into a hnite number of time steps with to = 0 < < ■ ■ ■ < = T 

and advance the solution forward in time. The transient term is typically 
discretized by an implicit scheme due to stability reasons which results in a 
(nonlinear) implicit residuum equation at each time step: 

R{y\y^-\y^-^...,u)=0 = (18) 

where y* ~ y{ti) approximates the solution at time ti. Depending on the 
order of the implicit time-stepping approximation, the residuum equation for 
a certain time step approximation y^ contains approximations to the solution 
at one or more previous time steps y^~^,y^~‘^ etc. Eor notational reasons, we 
choose the first order Backward Euler discretization and drop previous states 
except y''~^. Application of the proposed method to higher order schemes 
is straightforward. In case of the Backward Euler method, the residuum 
equations are 

R{y\ f-\u) := f ' - My\ u)=0 V* = 1,..., iV (19) 
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where fh represents a spatial discretization of the right hand side of the PDE. 

The set of nonlinear residuum equations can be solved one after another 
marching forward in time. Typically, iterative methods are used to converge 
to a pseudo-steady state solution at each time step: 


for i = 1,..., : 

iterate 2/fc+i = y" (20) 

with ?/q := where the iterator G* is designed such that the converged 
pseudo-steady states ?/* satisfy the residuum equations (19). We therefore 
assume, that G* is contractive with respect to ?/*, i.e. 

dG\y\y^-\u) 


dy^ 


<p*<l \/ i = 1,..., N 


( 21 ) 


for all points of interest, which ensures convergence of the above iterations. 
The converged pseudo-steady states are hxed points of G* such that y^' = 
holds if and only if R{y\y'^~^,u) = 0. 

In order to extend from simulation to single-step one-shot optimization, 
where one incorporates design updates already during the primal flow com¬ 
putation, the time-marching scheme (20) is modihed in such a way, that the 
residuum equations at each time step are solved inexactly. Instead, an outer 
loop is performed that updates the state at all time steps: 


iterate fc = 0,1... : 

VM=G‘(yl,vt+\,n) for i = (22) 

with y^ := y^ \/k G N. In contrast to ( [20| ), where hxed point iterations 
are performed at each time step to reach the converged states ?/* one after 
another, in the one-shot framework (22) a complete trajectory of the unsteady 
solution is updated within one iteration. Interpreting the time-dependent 
state variable as a discrete vector from the product space y = {y^,, y^) G 
YN .^y X ■■■ X Y we can write ( [^ in terms of an update function 

yk+i= (23) 

where H : x t/ —)■ Y^ performs the update formulas (22) and is dehned 

as 

( G\y\y\u) \ 


H{y,u) : = 


G^{y^G\y\y^,u),u) 

\G^{y^, G^-\y^-\ G^-\y^-^ ... ,G\y\y^,u),u)... ,u),u) J 


(24) 
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It is shown in [23], that the recursive iteration function H is contractive 
with respect to y, i.e. 


dH{y,u) 

dy 


<p<l 


(25) 


for all points of interest with p := maxjp*. This ensures convergence of the 
modihed time-marching scheme (22) to the hxed point y = H{y,u) where 
= G^{y\y^~^,u) holds for all i = By construction of H, the 

hxed point satishes the residuum equations (19) and is therefore a numerical 
solution of the unsteady PDE. 


3.3. Single-step one-shot optimization steps 

The proposed one-shot method aims at solving the following discrete op¬ 
timization problem 


min Jtv(?/,«) s.t. y = H{y,u) (26) 

y,u 


where Jat is an approximation of the time averaging objective function as for 
example 


JN{y,u) 


^^AtiJ{y\u) ^ 


- / J{y{t),u)dt 


(27) 


with Ati \= ti — ti-i. We dehne the corresponding Lagrangian function 


L{y, y, u) := J 7 v(?/, u) + {H{y, u) - yfy 


N 


where the adjoint variable is an element from the product space y E iY 
The necessary optimality conditions for the discrete optimization problem 
(26) yield the state, the adjoint and the design equation, each of which is 
analyzed in the sequel: 


1. State equation-. Differentiation of L with respect to the adjoint variables 
yields the unsteady residuum equations at each time step. As shown in 
the previous subsection, these equations can be solved simultaneously 
in one iteration that updates the state on the entire time domain: 


iterate A; = 0,1... : 

2/fc+i = H{yi^,u). 


( 28 ) 
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2. Adjoint equation: From VyiL = 0 for all i = 1,..., iV, the adjoint 
equation is derived: 


y = VyJN{y,u) + — 


dH\ 


- 


dy J 


y 


Since the transpose of the constraint Jacobian is an upper triangular 
matrix, the adjoint flow is backwards in time while := 0 is imposed 
at the last time step instead of initial conditions. In classical reduced 
space methods, the adjoint equations are solved iteratively one after 
another backwards in time. However, since we want to integrate the 
iteration into an one-shot optimization process, we rather solve the 
equations simultaneously for all time steps: 


iterate k = 0,1 


^ j , ^ , f dH{y,u)Y - 

Vk+I = ^yJNiy, U) + I —— 1 yf,. 


(29) 


The adjoint iteration (29) can be generated automatically by applying 


AD to the primal iteration and evaluating the objective function. 

3. Design equation: The design equation refers to stationarity of the La- 
grangian with respect to design changes and reads 

0 = VuJN{y,u) + i—j y. 

It is solved iteratively by 
iterate k = 0,1... : 

Uk+i =Uk- {VuJN{y, Uk) + duH{y, Ukfy) (30) 


with a preconditioning matrix Bk- 

In the single-step one-shot method, the set of all three equations is solved 
simultaneously in one single coupled iteration: 


'yk+1 


H{yk,Uk) 

yk+i 

= 

^yJN{yk^ Uk) + dyH{yf., Uk)'^yk 

Uk+1_ 


Uk - Bjj^ (v„ Jjv(?/fc, Uk) + duH{y^, Uk)'^yk)_ 


( 31 ) 


12 








As in the steady case, the matrix Bj. ensures convergence of the single-step 
one-shot method. In order to approximate the Hessian of the corresponding 
augmented Lagrangian (13) with respect to design changes, Bk is updated 
using BFGS-updates on the gradient VuL°‘ in each iteration. However, in 
contrast to single-step one-shot optimization with steady state PDEs, the 
update of the primal and adjoint iteration now each involve a loop over the 
entire time domain. Only this makes it possible to compute an approximation 
of the reduced gradient which is used in the design update. 

Enhancing a standard simulation code for single-step one-shot optimiza¬ 
tion involves only minor changes to the time-marching scheme. Since the 
inner iterator G* is used in a black box manner, the stability and robustness 
properties of the CFD solver are preserved within the new scheme. 


4. Single-step one-shot optimization with unsteady RANS 

The proposed one-shot method is applied to an optimal active flow con¬ 
trol problem of unsteady flow around a 2D cylinder. The flow is governed 
by the unsteady incompressible Reynolds-averaged Navier-Stokes equations 
(unsteady RANS) which are solved on a block-structured grid with 12640 
control volumes. 15 actuation slits are installed on the surface of the cylin¬ 
der where pulsed actuation is applied according to 

:= sin {2nft) — v} for / = 1,..., 15 (32) 


while the frequency / is hxed. The amplitudes v} at each slit are the design 
variables. The optimization objective is to hnd optimal actuation parameters 
u = (m^, ... ,M^®) that reduce vorticity downstream the cylinder. 

The optimization problem under consideration is given by 


dv 

m 


min 

y,u 


y / Cd{y{t),u)dt 


subject to 


uAv + {v ■ V)v + Vp = g 

divn = 0 
v{x, t) = (a^,..., 
v{x, t) = 0 
v{x, 0) = 


in D X (0, T] 

in D X (0, T] 
on F X (0, T] 
ondn\T X (0,T] 
in D 


(33) 


(34) 
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where Cd(t) denotes the drag coefficient around the cylinder. The state 
function y contains the velocity v{x, t) and the pressure p{x, t) in the domain 
n for a given force held g and viscosity z/ > 0. The actuation is applied on 
15 slits T distributed around the cylinder surface. 

The governing how equations are solved with the industry standard CFD 
simulation code ELAN [SU] . ELAN is a second order hnite volume code that 
approximates the transient term with the implicit BDF-2 scheme. The result¬ 
ing implicit equations are solved one after another using a SIMPLE scheme 
variant [31]. The SIMPLE scheme is a widely used numerical method for 
solving pressure-linked equations. In that scheme, pressure correction steps 
to the velocities are performed iteratively until a pseudo-steady state at that 
time step is reached. Notice, that all numerical ehort that makes the simu¬ 
lation code a stable and robust CFD solver is contained inside the pressure- 
correction iterations while the outer loop shifts the computed solution in 
time according to the transient approximation. We therefore identify the in¬ 
ner hxed point iterator G* as one step of the SIMPLE iteration at time step 
ti- The fixed point iterator H, that is used in the proposed one-shot method, 
performs a loop over all time steps applying one SIMPLE step at each time 
step. 

The AD tool Tapenade [32| is applied to the iteration function H as 
well as evaluating the discretized objective function. Its reverse mode auto¬ 
matically generates an iterative procedure for solving the adjoint equation 
and computing the reduced gradient. The design space preconditioner Bk 
is approximated using BFGS updates based on the reduced gradient, i.e. 
a = P = 0. The one-shot iteration stops, when a certain tolerance on the 
reduced gradient is reached || Jj -l- H'^vW < e, while in the present test case 
e = 0.001 was chosen. 

Figure [T] plots the residuals, the norm of the reduced gradient as well as 
the objective function during the single-step one-shot optimization. An av¬ 
erage drag reduction of about 30% is achieved with the optimization. Primal 
and adjoint residuals are reduced simultaneously with the reduced gradient 
indicating the successful application of the single-step one-shot method to the 
unsteady PDE-constrained optimization problem. The initial and optimized 
actuation is visualized in Figure]^ 

In this test case, the iteration number needed for convergence of the 
single-step one-shot method was observed to be only 3 times the number 
of iterations needed for a pure simulation with the modihed time-marching 
scheme with hxed design. This proves the typical bounded retardation of the 
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single-step one-shot method as it was observed numerically for the steady case 
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Figure 1: Optimization history of one-shot iterations solving unsteady incompressible 
Navier-Stokes equations (Re=100) 
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Figure 2: Initial and optimized amplitudes of the actuation. 


5. Improving primal convergence 

Since the bounded retardation property of the single-step one-shot method 
ensures that the cost for an optimization is only a small factor of the cost of 
a pure simulation, it is crucial to further investigate the performance of the 
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simulating fixed point solver. We therefore focus on an improvement of the 
primal iteration by applying adaptive time scales to the state variable. 

The hxed point iterations, that solve the unsteady PDE in the proposed 
one-shot method, perform one step of the inner update function G* at each 
time step. Due to the inexact approximation of states at previous time 
steps, the update at a certain time is contaminated by this error. This is 
also reflected in the lower triangular structure of the Jacobian of H. The 
errors are propagated through the entire time domain and accumulate until 
the last time step is reached. Thus, the number of iterations needed to 
reduce the residuals increases as the time domain is enlarged. It was observed 
numerically, that the dominating contribution to the error occurs in the 
direction of time while errors in the amplitude of the inexact trajectories 
are rather marginal. The computed trajectories exhibit a numerical time 
dilation compared to the hnal solution approximation. 

To reduce the time dilation and improve the primal convergence, we apply 
an adaptive time scaling approach as introduced in [23] • After each primal 
update, we assign a trajectory to a scaled time t such that the new trajectory 
is in phase with the physical solution. More precisely, we dehne the new 
trajectory as 

f-.= y{U) V* = l,...,iV (35) 


where ti is chosen such that the residual equation is minimized: 

.i—1 


mm 

u 


y - y 


ti—1 


- fh{y\u] 


Vi = 1,...,W 


(36) 


The global minimizer of (36) is given by 

{yi -yi-\fh{y\u)) 


ii — ti-l 


\\fk{v‘,u)\\\ 


(37) 


In this adaptive time scaling approach, we eliminate the error component 
in the direction of time in such a way, that the numerical time dilation 
vanishes. The convergence of the primal iteration is guaranteed with 


Wvk+i-yJl < Wvk+i-y. 

= \\Hiyk,u) ■ 


y. 


k^oo 


for any design u eU and the hxed point = H{y^,u). 


(38) 

(39) 
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6. Numerical results 

In this section, we investigate the performance of the fixed point iterator 
H and the effect of the adaptive time scaling approach. As a hrst test case, 
we consider a nonlinear ordinary differential eqnation (ODE), namely the 
Van-der-Pol oscillator. Since any unsteady PDE transforms into a system 
of ODEs after spatial discretization, the Van-der-Pol oscillator is used as a 
simple model problem. In order to take a step closer to the Navier-Stokes 
equations, we choose the one-dimensional linear advection-diffusion equation 
with periodic boundary conditions as a second test case. 


6.1. Van-der-Pol equation 

The Van-der-Pol oscillator is a nonlinear oscillator where a damping factor 
M > 0 controls the magnitude of the nonlinear term. It can be written as a 
system of two hrst order ODEs. With y = {x,v)'^ the Van-der-Pol oscillator 
reads 

\-x(t) + u{l - x(ty)v(t) 

a;(0)\ ^ 

^( 0 ); \Vo) 


Vt e (0, T] 


(40) 


where x and v denote the position and the velocity of the oscillator, respec¬ 
tively. 

Since we want to resemble the situation where the user is provided with 
an implicit time-stepping simulation tool, we approximate the transient term 
with the implicit Backward Euler method. The implicit equations are then 
solved one after another using an iterative Quasi-Newton method at each 
time step. According to Section the contractive function H is set up 
to converge the primal variable simultaneously for all time steps, while G* 
represents one step of the Quasi-Newton solver. To remove the numerical 
time dilation and improve the convergence behavior, the iteration function 
H is enriched with the adaptive time scaling approach according to (37). 

Figure visualizes the effect of the time scaling approach on the x- 
component of 16 intermediate trajectories during the primal computation. 
On the left, the numerical time dilation can be observed from the spanning 
bandwidth of intermediate trajectories. In contrast, the time scaled trajec¬ 
tories (right) are all in phase with the physical solution marked by triangles. 
The corresponding residuals are plotted in Figure]^ where the diagonal lines 
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Figure 3: cc-component of the Van-der-Pol oscillator for 16 different iterations without (left) 
and with time scaling approach (right). The physical solution is marked by triangles. 




Figure 4: Residuals of the Van-der-Pol oscillator for 16 different iterations without (left) 
and with time scaling approach (right). 


on the left indicate the dependency of the convergence on the time step 
number. When applying the time scaling approach (right), the residuals 
drop constantly over the entire time domain. This is also reflected in Figure 
1^ where the number of iterations that are needed for convergence are plotted 
for increasing numbers of time steps. Adapting the time scales in each iter¬ 
ation dramatically increases the performance of the primal iteration in this 
test case. 
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Figure 5: Number of iterations needed for solving Van-der-Pol oscillator with and without 
time scaling approach. 


6.2. Advection-diffusion equation 

As a second test case we investigate the advection-diffusion equation with 
periodic boundary conditions: 


dtyit, x) + ad^yit, x) - yd^^yit, x) 

y{f),x) 

2 /(^, 0 ) 


= 0 Vx e (0,l),t e (0,T] 

= h{x) Vx e [0,1] (41) 

= t/(t, 1) Vt>0 


with a = 1, /i = 10“® and initial condition h{x) = sin(27rx). We use the 
Backward Euler discretization in time with initial At = 0.01 and central hnite 
differences in space to approximate the derivative terms. A Quasi-Newton 
fixed point solver is applied at each time step to solve the implicit equations 
and resemble the scenario of a general implicit time-marching simulation tool. 
The iteration function H then performs one update of the Quasi-Newton 
solver at each time step according to section 

Figure 1^ plots 16 intermediate state trajectories at x = 0.5 while solving 
(41) with the hxed point iterator H. The numerical time dilation is obvious 
on the left and enlarges with time. Applying the time scaling approach 
(right) removes the time dilation such that all intermediate trajectories are 
in phase with the physical solution. The corresponding residuals are plotted 
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Figure 6: State of the advection-diffusion equation at a; = 0.5 for 16 different iterations 
without (left) and with (right) time scaling approach. Physical solution is bold. 




Figure 7: Residuals of the advection-diffusion equation for 16 different iterations without 
(left) and with time scaling approach (right). 


in Figure where the effect of the time scaling approach is visible on the 
right: The residuals drop with two orders of magnitude for all time steps 
after the first time scaling. Nevertheless, the dependency of the convergence 
on the number of time steps can still be observed from the diagonal structure 
of the residuals. 

In this test case, applying the time scaling approach still yields an im¬ 
provement for the primal iteration, but the independence of the iteration 
number for convergence from the number of time steps, as it was observed 
in the previous test case, could not be recovered for advection-driven flow. 
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7. Conclusion 

In the single-step one-shot method, the necessary optimality conditions 
for PDE-constrained optimization are solved simultaneously in the full space 
for the primal, the adjoint and the design equations. It is especially tailored 
for problems where the user is provided with a simulation tool that is to be 
used in a more or less black box fashion. Simulation tools for steady state 
PDEs often apply an explicit pseudo-time-stepping scheme which iterates in 
time until a steady state is reached. In the single-step one-shot optimization, 
these hxed point iterations are enriched by an adjoint and a design update 
step such that feasibility and optimality is reached simultaneously. The cost 
of an optimization with the single-step one-shot method has proven to be 
only a small multiple of one pure simulation. 

However, if the PDE is unsteady, setting up an appropriate fixed point 
solver is non-trivial since common schemes often apply an implicit time¬ 
marching method and solve the residual equations one after another forward 
in time. It has been shown in this paper, that these time-marching schemes 
can be modihed to £t into the proposed one-shot optimization framework by 
reducing the number of inner iterations for solving the implicit equations. 
In the resulting approach, the entire time trajectory of the unsteady PDE is 
updated within one iteration. Applying AD to the modified time-marching 
scheme as well as evaluating the discrete objective function automatically 
generates a consistent discrete adjoint iteration. Augmenting the modified 
primal iteration with an adjoint update and a preconditioned reduced gradi¬ 
ent step for the design yields the single-step one-shot optimization method 
for unsteady PDE-constrained optimization. The method has been applied 
to an optimal active flow control problem governed by the unsteady RANS 
equations. 

The modified time-marching scheme has been further improved apply¬ 
ing adaptive time scales. In that approach, the intermediate trajectories are 
shifted in time such that they are in phase with the physical solution. Eor 
a general model problem, the time scaling approach yields an improved con¬ 
vergence behavior that is independent of the time domain resolution. An 
application to unsteady flow driven by advection showed that the time scal¬ 
ing approach still yields an improvement for the convergence speed while the 
dependence on the number of time steps remains. 
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